Evaluación de la pérdida de calidad del audio usando el codec OPUS

Juan M. Fonseca-Solís · Agosto 2020 · 5 min read


Resumen

Los estetoscopios digitales son herramientas que le permiten a los médicos realizar remotamente un seguimiento de sus pacientes por problemas de salud relacionados, por ejemplo, con el corazón y los pulmones. En muchos de los casos, estos estetoscopios digitales utilizan los mismos codecs que son empleados en sistemas de teleconferencia como Zoom, MS Teams y Skype. El éxito de estos codecs está en que son capaces de transmitir voz y música de alta calidad empleando bajas tasas de bits; para ello se valen de que el oído humano es incapaz de escuchar ciertas frecuencias y por eso las eliminan. Aunque esto resulte excelente para cuestiones de recreación y comunicación, no es tan bueno para los algoritmos de reconocimiento de patrones, que necesitan la información completa para sus algoritmos. En este ipython notebook nos dimos a la tarea de realizar pruebas al códec OPUS para determinar si la pérdida de la calidad es realmente significativa. Los resultados mostraron que, para un algoritmos de reconocimiento de patrones robusto y aplicando un filtrado pasabajas y una tasa de bits de 8kbps, el algoritmo OPUS parace ser apropiado para transmitir sonidos producidos por el corazón y los pulmones sin necesidad de transmitir datos crudos.

Códec OPUS

En comparación a otros algoritmos como el MP3, AAC y Vorbis, el códec OPUS ofrece una calidad de audio similar pero a una menor tasa de bits y es capaz de cubrir casi todo el rango de transmisión (ver figura 1). Esto significa que se puede ofrecer una calidad de audio mejor al mismo costo. Opus está compuesto por dos algoritmos llamados SILK (creado por Skype Limited para comprimir voz) y Constrained Energy Lapped Transform (CELT) (para comprimir música) [3, 4, 13]; así como un modo híbrido.

Figura 1. Tomada de https://opus-codec.org/static/comparison/quality.svg.

Psicoacústica

La mayoría de códecs de sonido utilizan la teoría psicoacústica para "engañar" al oído humano (más exáctamente, la cóclea) haciendo que perciba la misma calidad de audio a costa de eliminar las frecuencias inaudibles. Para ello, utilizan conceptos tales como las escalas Bark y ERB, bandas críticas, enmascaramiento frecuencial, enmascaramiento temporal y predicción lineal (LPC, por sus siglas en Inglés). A continuación hacemos un resumen de estos conceptos:

  • Cóclea: órgano en forma cónica (ver figura 2) encargado de percibir las frecuencias y transmitirlas como impulsos eléctricos al cerebro. Se encuentra "arrollado" dentro del oído interno y recibe estímulos del oído externo e intermedio mediante la ventana oval (esta se conecta con los huesecillos conectados al tímpano).

Figura 2. Tomada de https://en.wikipedia.org/wiki/Auditory_system#/media/File:Anatomy_of_the_Human_Ear.svg.

  • Escalas Bark y ERB: modelos matemáticos utilizados para representar la resolución frecuencial de la cóclea. Se basan en el hecho de que el oído humano realiza un mejor trabajo reconociendo tonos graves que tonos agudos (distribución logarítmica). La escala Equivalent Rectangular Bandwidth (ERB) es una versión mejorada de la escala Bark y su ecuación es $\text{ERB}(f)=21.4 \log_{10}{(0.00437f+1)}$, mostrada a continuación [9].
In [1]:
%pylab inline
import numpy as np

def erb(F):
    return 21.4 * np.log10(0.00437*F+1)

F = np.linspace(0,20000,1000) # 0 a 20k Hz es el rango de audición humana
E = erb(F)
pylab.plot(F/1e3,E)
pylab.ylabel('Frequency (kHz)')
pylab.xlabel('ERB (kHz)')
pylab.show()
Populating the interactive namespace from numpy and matplotlib
  • Bandas críticas: segmentos de la cóclea, de diferente longitud y localización, en los cuales el oído humano no es capaz de diferenciar dos tonos lo suficientemente cercanos. Dos ejemplos de bandas críticas se ubica en el rango $[555, 677]$ kHz con frecuencia central de 511 Hz y $[8928, 10353]$ Hz con frecuencia central 9589 Hz [10].
In [31]:
E = np.linspace(0,erb(20000),32) # 32 puntos equidistantes
F = (np.power(10,E/21.4)-1)/0.00437 # invertimos el logaritmo de la escala ERB
print(F.astype(int))
pylab.stem(F,np.ones(len(F)))
pylab.xlabel('Frecuencia central de cada subfiltro (Hz)')
pylab.show()
[    0    35    76   124   179   242   315   400   498   611   742   893
  1068  1270  1503  1772  2083  2443  2859  3339  3894  4536  5277  6134
  7123  8267  9589 11116 12881 14920 17276 20000]

Diseño de los casos de prueba

Basándonos en los conceptos de psicoacústica explicados arriba, y considerando que el rango de frecuencias para los sonidos producidos por los pulmones es de $[60,1000]$ Hz y del corazón es de $[20, 500]$ Hz, podemos definir una serie de casos de prueba para determinar si la pérdida de información frecuencial es significativa [14, 15].

  • Caso de prueba 1: el contenido de la magnitud espectral presente en un barrido de frecuencias de $[20,1000]$ debe mantenerse, sin que se agreguen o eliminen frecuencias.
  • Caso de prueba 2: para dos tonos con frecuencia 590 kHz y 610kHz (ubicados en la misma banda crítica), solo el de 610kHz deberá escucharse si la energía del de 590 kHz está por debajo de cierto umbral (enmascariento frecuencial).
  • Caso de prueba 3: para una señal real y cruda captada por un estoscopio digital, determinar si hay afectación al ser comprimida y descomprimida.

Para probar el caso 1 se puede utilizar un barrido de frecuencias que siga la siguiente ecuación [6] $$ F(t) = \Big(\frac{F_1-F_0}{T}\Big)t + F_0, $$

El rango de frecuencias a usar es de $[20, 20000]$ Hz (el rango de audición humana). El barrido puede ser construido usando un senoidal de la forma $x[t] = A \sin(2\pi Ft)$.

In [17]:
%pylab inline
from scipy.io import wavfile
from IPython.display import Audio
import numpy as np

def plotSpecgram(x,Fs,fMax=None):
    # spectrogram
    fig, ax = pylab.subplots(nrows=1)
    ax.specgram(x, NFFT=1024, Fs=Fs, noverlap=900)
    pylab.xlabel('Tiempo (s)')
    pylab.ylabel('Frecuencia (Hz)')
    if fMax != None:
        pylab.ylim([0,fMax])
    pylab.show()
    
def plotFFT(x,Fs,xlim=None):
    pylab.figure()
    N2=int(len(x)/2)
    F = np.linspace(0,Fs/2,N2)/1e3
    X = np.sqrt(np.abs(np.fft.fft(x)[0:N2]))
    pylab.plot(F, X)
    X[0] = 0 # remove DC value
    pylab.xlabel('Frecuencia (kHz)')
    pylab.ylabel('$\sqrt{|S(f)|}$')
    if xlim != None:
        pylab.xlim(xlim)
    pylab.show()

rango = [20.0, 20000.0] # el rango promedio de audicion humano en Hz
Fs = 44.1 * 1e3 # la tasa de muestreo de los equipos comerciales
T = 1.0 # segundos (t1-t0)

N = int(T*Fs)
n = np.arange(0,N)
F0 = (rango[1]-rango[0])*n/N + rango[0]

x = np.sin(np.pi*F0/Fs*n) # f=F0/Fs: frecuencia discreta

plotSpecgram(x,Fs)
plotFFT(x,Fs)

wavfile.write('/tmp/barrido.wav',int(Fs),x.astype(np.float32))
Audio(x, rate=Fs) 
Populating the interactive namespace from numpy and matplotlib
Out[17]:

El espectrograma muestra que el barrido va de los 20 a los 20000 Hz como se deseaba. Asimismo, la magnitud espectral muestra una recta horizontal casi en todo el rango, indicando que todas las frecuencias están presentes como se supone que tiene que pasar.

Para el caso de prueba 2, podemos construir una señal de dos tonos como sigue:

In [40]:
y = np.sin(np.pi*590/(Fs/2)*n) + np.sin(np.pi*610/(Fs/2)*n) 
plotSpecgram(y,Fs,fMax=2000)
plotFFT(y,Fs,xlim=[500/1e3,700/1e3])

wavfile.write('/tmp/enmascaramiento.wav',int(Fs),x.astype(np.float32))
Audio(y, rate=Fs)
Out[40]:

Ejecución de los casos de prueba

Caso de prueba 1

Procedemos a codificar y decodificar el barrido de frecuencias usando OPUS para una tasa de bits de 8 Kbps (la misma empleada en una videollamada entre dos personas) [12].

In [34]:
def readPlayVisualizeFile(inputFile,fMax=None):
    fs, x = wavfile.read(inputFile)
    y = np.array(x)/max(x)
    if None==fMax:
        plotSpecgram(x,fs)
    else:
        plotSpecgram(x,fs,fMax)
    return fs, x

#!sudo apt-get install ffmpeg
#!sudo apt-get install opus-tools
!ffmpeg -loglevel error -y -i /tmp/barrido.wav -qscale 0 /tmp/wavRaw.wav # opus-tools
!opusenc --quiet --bitrate 8 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav

Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)

Audio(x, rate=Fs)
Skipping chunk of type "LIST", length 26
Out[34]:

El espectrograma muestra que partir de los casi 5 kHz el barrido de frecuencia se pliega, en un efecto llamado aliasing, lo que produce que veamos frecuencias que en realidad no existen, esto sí que podría afectar la calidad del clasificador. Para resolverlo se podría emplear un filtro pasabajas con frecuencia de corte 5 kHz antes de realizar la compresión. Por su parte, la magnitud espectral muestra que a partir de los 5 kHz las frecuencias altas son atenuadas. Probemos ahora con 128 kbps (fullband stereo) a ver si la calidad mejora. Sin embargo, aunque se introdujeron artefactos que antes no estaban, el rango $[20,1000]$ Hz sigue estando casi intacto, lo cual indica que para un algoritmo de reconocimiento de patrones lo suficientemente selectivo, la corrupción en la señal no debería afecta (hay que considerar si existen armónicas para asegurar esto completamente).

In [35]:
!ffmpeg -loglevel error -y -i /tmp/barrido.wav -qscale 0 /tmp/wavRaw.wav # opus-tools
!opusenc --quiet --bitrate 128 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav

Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)

Audio(x, rate=Fs)
Skipping chunk of type "LIST", length 26
Out[35]:

El espectrograma sigue mostrando que hay cierto aliasing pero con bastante menor energía, lo cual favorable; y la magnitud espectral muestra una recta horizontal similar a la que vimos con la señal cruda. Por lo que podemos decir que con full-band ofrece una menor distorsión. En este punto valdría la pena preguntarse cuál de los tres métodos de OPUS realizó la compresión y si los resultados variarían en caso de que se fuerce algún método en particular (SILK, CELT o el híbrido).

Caso de prueba 2

Lorem ipsum...

In [39]:
!ffmpeg -loglevel error -y -i /tmp/enmascaramiento.wav -qscale 0 /tmp/wavRaw.wav # opus-tools
!opusenc --quiet --bitrate 8 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav

Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)

Audio(x, rate=Fs)
Skipping chunk of type "LIST", length 26
Out[39]:

Caso de prueba 3

Nota: a partir de este momento se recomienda usar audifonos para apreciar mejor la calidad el audio.

Procedemos ahora a codificar y decodificar una grabación real y cruda usando OPUS para una tasa de bits de 8 Kbps (la misma empleada en una videollamada entre dos personas) [12]. Procedimos a tomar tomar una grabación del sitio https://www.kaggle.com/vbookshelf/respiratory-sound-database, la cual corresponde a una grabación de tomada de un individuo .

In [32]:
Fs, x = readPlayVisualizeFile('./wav/107_3p2_Tc_mc_AKGC417L_2.wav')
plotFFT(x,Fs)

Audio(x, rate=Fs)
Out[32]:

Observamos que hay energía alta en el rango $[0,1.5]$ y otra menor en $[7.5,9.5]$ kHz. La primera parece corresponder a los sonidos naturales del cuerpo, la segundo para más bien una frecuencia parásita que podría haberse colado al momento de realizar la grabación (el cuerpo humano no produce tonos constantes y puros como los que se observan después de los 7500 Hz).

In [33]:
# OPUS
!ffmpeg -loglevel error -y -i ./wav/107_3p2_Tc_mc_AKGC417L_2.wav -qscale 0 /tmp/wavRaw.wav # same quality
!opusenc --quiet --bitrate 8 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav

Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)

Audio(x, rate=Fs)
Skipping chunk of type "LIST", length 26
Out[33]:

El espectrograma muestra que la energía un poco antes de los 5 kHz fue atenuada, pero ya descubrimos que esto no importa tanto pues se sale del rango objetivo en los $[20,1000]$ Hz. Lo más importante parece ser que las frecuencias están inalterada en el rango deseado. La magnitud espectral muestra un espectro limpio.

4 Conclusiones

  • La calidad del OPUS parece aceptable para la tasa de bits usada en las llamadas 1:1 (8kps) cuando los sonidos a analizar corresponden a los del corazón y los pulmones.
  • Antes de emplear OPUS con una tasa de bits de 8kps, es necesario aplicar un filtro pasabajas con frecuencia de corte 5 kHz para evitar el efecto del aliasing y confundir al algoritmo de reconocimiento de patrones.
  • Emplear el OPUS bajo estas dos condiciones explicadas en los puntos anteriores parece ser suficiente y no es necesario enviar los datos sin comprimir.
  • Como trabajo futuro, se podría investigar un método para cuantificar la distorsión apreciada en los casos de prueba propuestos y contemplar más casos de prueba (por ejemplo, para el enmascaramiento temporal y evaluar una colección de sonidos).

Referencias

  1. Colaboradores de Wikipedia. Códec de audio [en línea]. Wikipedia, La enciclopedia libre, 2020 [fecha de consulta: 4 de noviembre del 2020]. Disponible en https://es.wikipedia.org/w/index.php?title=C%C3%B3dec_de_audio&oldid=129362034.
  2. Hong Kong Polytechnic University. Department of Electronic and Information Engineering. Perceptual Coding and MP3. Disponible en http://www.eie.polyu.edu.hk/~enyhchan/DAP-lab-PModel-v3.pdf.
  3. López Monfort José Javier. Opus codec | 22/23 | UPV. Universitat Politècnica de València - UPV. Disponible en https://youtu.be/2-yv1bCDL94.
  4. Wikipedia contributors. (2020, October 11). SILK. In Wikipedia, The Free Encyclopedia. Retrieved 15:18, November 8, 2020, from https://en.wikipedia.org/w/index.php?title=SILK&oldid=982904860
  5. Wikipedia contributors. (2020, November 2). Linear predictive coding. In Wikipedia, The Free Encyclopedia. Retrieved 15:22, November 8, 2020, from https://en.wikipedia.org/w/index.php?title=Linear_predictive_coding&oldid=986664067
  6. EPFL. Room impulse responses. Ultima vez consultado el 27 Dec 2020 en: https://nbviewer.jupyter.org/github/LCAV/SignalsOfTheDay/blob/master/Room_Acoustics/Room%20Impulse%20Response.ipynb
  7. López Monfort José Javier. Enmascaramiento Temporal | 4/23 | UPV. Disponible en https://youtu.be/7NQvDoZMBM8 (consultado por última vez el 3 de Enero de 2021).
  8. International Telecomunication Union (ITU-T). P.800.1, MOS.
  9. Julius O. Smith III, Jonathan S. Abel. Equivalent Rectangular Bandwidth. Center for Computer Research in Music and Acoustics (CCRMA), Stanford University. Disponible en https://ccrma.stanford.edu/~jos/bbt/Equivalent_Rectangular_Bandwidth.html (consultado por última vez el 5 de Enero de 2021).
  10. LCAV. Python MP3 encoder. Disponible en https://github.com/LCAV/MP3Lab/blob/master/mp3python.ipynb (consultado por última vez el 6 de Enero de 2021).
  11. XIPH ORG. OpusFAQ. Disponible en https://wiki.xiph.org/OpusFAQ#Why_not_keep_the_SILK_and_CELT_codecs_separate.3F (consultado por última vez el 6 de Enero de 2021).
  12. Tyler Abbott. How Much Data Does a Zoom Meeting Use?. Reviews ORG. Disponible en https://www.reviews.org/internet-service/how-much-data-does-zoom-use/ (consultado por última vez el 16 de Enero de 2021).
  13. Xiph.Org Foundation. The Opus Codec. 135th AES Convention2013 October 17–20 New York, USA. Disponible en https://arxiv.org/pdf/1602.04845.pdf (consultado por última vez el 21 de Enero de 2021).
  14. GROSS, V., DITTMAR, A., PENZEL, T., SCHÜTTLER, F., & von WICHERT, P. (2000). The Relationship between Normal Lung Sounds, Age, and Gender. American Journal of Respiratory and Critical Care Medicine, 162(3), 905–909. https://doi.org/10.1164/ajrccm.162.3.9905104 (consultado por última vez el 23 de Enero de 2021).
  15. McGee, S. (2018). Chapter 39 - Auscultation of the Heart: General Principles. In S. McGee (Ed.), Evidence-Based Physical Diagnosis (Fourth Edition) (Fourth Edition, pp. 327-332.e1). Elsevier. https://doi.org/https://doi.org/10.1016/B978-0-323-39276-1.00039-1 (consultado por última vez el 23 de Enero de 2021).

Licencia de Creative Commons
Este obra está bajo una licencia de Creative Commons Reconocimiento-NoComercial-SinObraDerivada 4.0 Internacional. El sitio juanfonsecasolis.github.io es un blog dedicado a la investigación independiente en temas relacionados al procesamiento digital de señales. Para reutilizar este artículo y citar las fuente por favor utilice el siguiente Bibtex:

@online{Fonseca2020,
  author = {Juan M. Fonseca-Solís},
  title = { Pruebas por pares o pairwise testing},
  year = 2020,
  url = {https://juanfonsecasolis.github.io/blog/JFonseca.pairwisetesting.html},
  urldate = {}
}